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ABSTRACT 

Binary systems emit gravitational waves in a well-known pattern; for binaries in 
circular orbits, the emitted radiation has a frequency that is twice the orbital 
frequency. Systems in eccentric orbits, however, emit gravitational radiation in the 
higher harmonics too. In this paper, we are concerned with the stochastic background 
of gravitational waves generated by double neutron star systems of cosmological origin 
in eccentric orbits. We consider in particular the long-lived systems, that is, those 
binaries for which the time to coalescence is longer than the Hubble time lOGyr). 
Thus, we consider double neutron stars with orbital frequencies ranging from 10“® to 
2 X 10“®Hz. Although in the literature some papers consider the spectra generated 
by eccentric binaries, there is still space for alternative approaches for the calculation 
of the backgrounds. In this paper, we use a method that consists in summing the 
spectra that would be generated by each harmonic separately in order to obtain 
the total background. This method allows us to clearly obtain the influence of each 
harmonic on the spectra. In addition, we consider different distribution functions for 
the eccentricities in order to investigate their effects on the background of gravitational 
waves generated. At last, we briefly discuss the detectability of this background by 
space-based gravitational wave antennas and pulsar timing arrays. 
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1 INTRODUCTION 

It is well known that binary systems emit gravitational 
waves (GWs) and, for circular orbits, the frequency of the 
emitted radiation is twice the orbital frequency. On the other 
hand, for eccentric orbits, the pattern of emission has a 
more complex form: the emitted signals are also formed by 
the higher harmonics of the orbital frequency. In fact, the 
higher the orbital eccentricity is, the greater the contribution 
of the higher harmonics to the emitted radiation (se e e.g., 
IPeters fc Mathe^ 1 196,31 : IPouglass fc Braginskv^fI979^ . 

Since ma ny binaries can have non-negligible orbital 
eccentricities l|Kowaska et al.l l2012l l. it is interesting to 
investigate the behaviour of the spectra when eccentricity is 
taken into account. Particularly, this paper is addressed to 
the stochastic background of GWs generated by long-lived 
double neutron star (DNS) systems of cosmological origin 
in eccentric orbits. We consider that a system is long-lived 
if its lifetime is longer tha n the Hubble time, namely ~ 10 
Gyr (iKowaska et al. l2012h . Quantitatively, this assumption 
implies that the binaries have orbital frequencies in the range 
10“® — 2 X 10“®Hz, as explained in Section[5] We should bear 
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in mind that in this case we can neglect the time evolution 
of the orbital parameters of the systems. 

It is worth noting that the inclusion of eccentricity 
in the calculation of stochastic backgrounds generated by 
compact binaries is not new. In fact, one can find in the 
literature interesting papers on the stochasti c backgrounds 


2001 

). IKowaska et al.l ( 

2007 

i. In lignatiev et al.l 


and lEnoki fc Nagashima 


ll200lh , the authors study, for 
example, the GW background produced by coalescing 
neutron stars in t he Galaxy in eccentric orbits; in 
iKowaska et al.l (|2ni2l i the authors study the backgrounds 
generated by Population III binaries, giving special 
attention to the s pectra generated by l ong-liv ed and short¬ 
lived binaries; in lEnoki fc Nagashimal ([2003), the authors 
studied the background generated by supermassive black 
hole binaries. 

In general, the procedure to obtain the backgrounds 
generated by eccentric binaries includes, as a key 
component, the spectra of energy emitted by the source, 
which is written as dE/du, where E is the energy and u 
is, in our case, the orbital frequency. In the literature, the 
authors usually handle such a quantity following the steps: 
first, they sum the spectra dE/dVn corresponding to each 
harmonic, where n is the number of the harmonic, in order 
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to obtain the total spectrum dE/du\ subsequently, dE/du 
is included in the formula that determines the stochastic 
background. Here we use a different method: we calculate 
the backgrounds that would be generated by each harmonic 
of the emitted frequency and then we sum all these signals 
in order to obtain the total background. As we will see, this 
method allows us to clearly analyse the influence of each 
harmonic on the resulting backgrounds. 

In addition, we consider in our calculations different 
distribution functions for the orbital eccentricity e, and 
study the influence of such distributions on the resulting 
backgrounds. Particularly, we considered the case in which 
the distribution function of the eccentricities has the form 
of a Dirac delta function, that is, all the systems have the 
same eccentricity. We consider scenarios in which all systems 
have eccentricities of 0.25, 0.5, 0.75 and 0.9. Moreover, 
we consider the case of circular orbits (null eccentricity) 
and also another distribution that covers a wide range of 
eccentricities, namely, a uniform distribution ranging from 
e = 0 to e = 0.9. 

The paper is organized as follows: in Section [21 we 
derive the fundamental equations used to calculate the 
stochastic background of GWs generated by DNS systems 
of cosmological origin in eccentric orbits; in Section O we 
calculate the spectrum of energy for eccentric systems; in 
Section (4) we present and discuss the distribution functions 
of the eccentricities /(e) adopted; in Section [S] we discuss 
the maximum and minimum frequencies of the spectra and 
the population characteristics of the systems, such as the 
star formation rate density p'*( 2 ) and the mass fraction Adns 
of stars that is converted in DNSs; we present the results in 
Section and the conclusions and perspectives in Section (7] 


2 FUNDAMENTAL EQUATION 


For the calculation of the spectra generate d by each 
harmo nic, we follow the method describ ed in iRegimb^ 
II2OIID and iRegimbau fc Chauvineaul (l2007ll . who calculated 
the backgrounds in terms of the energy density parameter 
flgw, namely 


Hgw(r'o) 


- ttVoF{vo) 

pcC^ 


( 1 ) 


where pc = SHq /8ttG is the critical density of the universe. 
Ho = 70km s“^Mpc“^ is the Hubble parameter, vq is the 
observed frequency and F{vo) is the spectra l density of the 
flux, which is given by ([Ferrari et al. l ll999alH 


F{vo) = J f{vo,z)^dz, 


( 2 ) 


where 2 is the redshift and f{vo,z) is the energy flux per 
unit f requency generat ed by a unique source, which has the 
form (lRegimbaull201ll 'l 

f[uo,z) = —^^ . (3) 

47rr(2)2 dv v=vo{i+z) 

In equation (0, r( 2 ) is the comoving distance, v = 
vq{1-\-z) is the emitted frequency and dE/du is the spectrum 
of energy emitted by a binary system, calculated in the 
source frame. In addition, in equation ®, the event rate 
dR/dz is calculated by 


dR _ ^ p,{z)dV 
dz ^ 1 + 2 d 2 ’ 


( 4 ) 


where p*{z) is the star formation rate density (SFRD), Adns 
is the mass fraction of stars that is converted into DNSs, 
and dV/dz is given by 




rjzf 

■\/Hm (1 + 2)^ + Hk(l + 2)^ + Ha 


( 5 ) 


We use equation in a scenario in which the density 
parameters obey the relations Hm = Hdm + He and Hk + 
Hm + Ha = 1, with the subindices M, DM, B, k and 
A corresponding to matter, dark matter, baryonic matter, 
curvature and cosmological constant, respectively. Here, we 
consider a cosmological scenario in which Hm = 0.3, Hb = 
0.04, Hk = 0 and Ha = 0.7. 

In addition, the eccentricity is introduced in the 
calculation of the background by means of a distribution 
function /(e) in equation ([!]). So, with the above 
considerations and performing all the substitutions, 
equation 0 assumes the following form 


o t,. ^ 


/ 


3773 C2 ' 




pjz) 


dE 


no) + AE{z) dv 


U = Uo(l + 2) 


fie)dz, 


( 6 ) 


where the meaning of the limits of integration, Zf{vo) and 
Zi{vo), will be considered in Section[5l 


3 ENERGY EMITTED BY ECCENTRIC 
BINARIES 


In the previous section, we noticed that the spectrum of 
energy dE/dv is a key ingredient for the calculation of the 
backgrounds, since it provides information concerning the 
patterns of emission of the systems as a function of the 
orbital frequency. We can now obtain an expression for the 
energy spectrum for eccentric systems starting from 

dE dE dt dH 

dvn dt dfl dvn ’ 


where we included the index n in the notation, since 
equation 0 refers to the spectrum generated by each 
harmonic. In addition, H is the angular orbital frequency 
and dt/dil is given by 

dH dH da , , 

dt da dt ’ 


in which ’a’ represents the orbital semi mafor axis. In our 
case, dE/dt and da/dt are given by ([Peters &: Mathews! 
Il963ll 


dE 32 mfm|(mi + m 2 ) 


dt 


5 c® 


5(w,e); 


( 9 ) 


and 

da 

dt 


64 G^mim2(mi + m 2 ) 
5 c®a® 


FiA, 


( 10 ) 


where mi and m 2 are the masses of the components of the 
system. In addition, the functions g(n,e) and F{e) have the 
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form 

2 

Jn- 2 {ne) — 2eJ„-i(ne) H— Jn(ne)+ 
n 

+2eJn+i(ne) - Jn+ 2 ine)f + (1 - e^) [Jn- 2 {n£) + + 

-2£Jn{n£) + Jn+2ine)f + , ( 11 ) 

where Jn{ne) are the Bessel functions of order n, and 



F{e) = 


1 + 73eV24 + 37£'‘/96 

(1- £2)7/2 


( 12 ) 


Moreover, equations m and m are related to each 
other by 


F(£) = ^g(n,£), (13) 

n=2 

and dfl/da is obtained from Kepler’s third law. 

Now, recall that systems in eccentric orbits emit 
radiation at harmonics of the orbital frequency, which is 
related to the angular orbital frequency by = 27rr'orb. 
Thus, we can write the emitted frequency Vn of the harmonic 
n as a function of n and fl, namely 

TlO 

i^n = — for n ^ 2, (14) 

ZTT 

which enable us to calculate dQ/dVn in equation ([7|l. 

Now, performing all the necessary substitutions, 
equation 0 is transformed into 

dE _ 1 / 27rG \^/^ mim2 -i/3 g{n,e) , . 

du„ 3\nJ (mi+m2)i/® " T’(e) 

Note that for n = 2, E(£) = 1 and g{n,£) = Sn 2 , 
equation m c omes down to the ex pression valid for circular 
orbits (see e.g., [Ferrari et aLll200lh . 

In addition, g{n,£)/F[£) has the following meaning: it 
provides the fraction of energy emitted at the harmonic n 
by a system with orbital eccentricity £. 


4 THE DISTRIBUTION FUNCTIONS OF THE 
ECCENTRICITIES 


We are interested in studying the role of the distribution 
function of the orbital eccentricities on the stochastic 
background of GWs generated by DNS systems of 
cosmological origin. We consider distributions that, 
although could be considered unrealistic, can provide an 
insight into the influence of these functions on the spectra. 

First, we consider a scenario where all the DNSs have 
the same orbital eccentricity, such that 

f{£)=S{£-£'), (16) 


where 5{£ — S) is a Dirac delta function. Particularly, we 
consider the cases for e' = 0 (a population of circular DNSs), 
0.25, 0.5, 0.75 and 0.9. 

We also consider a distribution that covers a wide range 
of values of £. For this purpose, we use a uniform and 
normalized function that is given by 


/(£) = 


( 


(s^max S^min ) ^ for 


^ e ^ £n 


0 otherwise. 


(17) 


where Emin = 0 and Emax = 0.9 are the minimum and 
the maximum eccentricities, respectively. This distribution 
provides an interesting case to be compared with the 
distribution functions given in terms of the Dirac delta 
functions. 


5 THE STOCHASTIC BACKGROUND 


Before proceeding, we pay attention to the choice of the 
minimum and maximum values of the orbital frequencies 
used in the calculations of the backgrounds. First, as already 
mentioned, we consider in this paper long-lived DNSs, which 
lead us to the determination of a maximum value for the 
initial frequency. We mean by a ’long-lived system’, a binary 
whose lifetime is at least of the order of the Hubble time, 
that is, ~ lOGyr. 

In addition, recall that the greater the orbital 
eccentricity is, the shorter is the binary system lifetime. 
Therefore, we choose the maximum frequency such that 
the most eccentric systems considered in this paper have 
lifetimes Tufe equal to the Hubble time. The expressio n for 
Tiife, considering e close to 1, is given by (jPeterallQBdh 


Tlife 


768 Qq _ -2\7/2 
425 4,3 ’ 


(18) 


where fi = 64G®mim2(mi-|-m2)/5c®, ao is the initial orbital 
semimajor axis and eq is the initial eccentricity. Considering, 
for example, Tufe ~ lOGyr, mi = m 2 = 1.4M0 and £0 = 0.9 
in equation (1181) and using Kepler’s third law, we obtain a 
frequency of « 2 x 10“®Hz. 

On the other hand, according to lRosa d3 (l201l[) . there is 
a minimum frequency one must consider in the study of GWs 
generated by binary systems. He argues that the emitted 
radiation with frequencies lower than a given value should 
not be considered, since other mechanisms th at take away 
energy from the systems are more effective (see I Rosa d3l201ll 
for a detailed discussion.) _ _ 

Still following the discussion of iRosa 13 ll201lD , we 
should bear in mind two facts: it is difficult to define exactly 
the values of the minimum frequencies; on the other hand, 
such values are not important in practice, since the choice of 
a particular frequency will not have noticeable effects on the 
results. Gonsidering these arguments, we adopted the rather 
arbitrary value of 10“®Hz in our calculations. 

To proceed, the next step is to determine the SFRD 
Pt(z) and the value of Adns- 

Goncerning the SFRD, many different proposals are 
considered in the literature, although they do not differ from 
each other very si gnificantly. Here we adopt, as a fiducial 
one, that given bv ISpringel fc Hernauisd (l2003h . namely 






(19) 


where a = 3/5, fi = 14/15, Zm = 5.4 and pm ~ 
O.15M0yr“^Mpc“^ fixing the normalization. It is worth 
mentioning that the authors considered a AGDM cosmology, 
where the density parameters have the values given in 
Section [ 2 ] Moreover, they considered the Hubble parameter 
as Hq = 100 h km s“^ Mpc“^ with h = 0.7. 

_ The mass fraction Adns (see e.g.. iRegimbau fc Pache^ 

l2006l. for a detailed discussion) is given by 
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Ad ns — Pnsfp^ n 


( 20 ) 


where /3ns is the fraction of binary systems that survive 
to the second supernova event, fp provides the fraction 
of massive binaries (that is, binary systems where both 
components could generate a supernova event) formed inside 
the whole population of stars, and >I>ns is the mass fraction 
of neutron star progenitors. In our case, we assume that "hns 
is calculated by 


<hns = 




( 21 ) 


where is the Salpeter mass function llSalpeterjllOSSll . 

with (f){rn) = A = 0.17 and x = 1.35. Moreover, 

note that we are considering that the progenitors of the 

neutron stars have masse s in the range 8 — 25M |3. _ 

According to iRegimbau fc Pacheco! ll2006ll . 
the parameters /3ns and fp depend on different assumptions. 
First, in the evolutionary scenario of massive binaries the 
authors considered that none of the stars underwent recycle 
accretion. On the other hand, the velocity distribution of the 
natal kick affects the parameter /3ns, since an imparted kick 
could disrupt binaries or, less likely, prevent the disruption 
of the system. Particularly, the authors adopted a one¬ 
dimensional velocity dis persion of ss 80km s~^ . _ 

Still according to iReeimbau fc Pacheco! (l2006l ~). the 
estimates of fp depend on /3ns itself and on the ratio 
between single and double NS systems in the Galaxy. As 
a result, t here are uncertainties in th e values of these 
parameters. IRegimbau fc Pacheco! ll2006h then estimate the 
relative uncertainties in fp and /3ns. Numerically, they 
obtained o-fpifp « 0.5 and <Jp-^^/Pns « 0.75, leading to a 
relative uncertain of crAj£,,g/Adns ~ 0.9 in Adns. The authors 
stated that these uncertainties are solely formal, resulting 
from the simulations they performed, which depend on the 
evolutionary scenario for the progenitors. 

Concerning the range of masses of the neutron star 
progenitors, recall that the minimum and maximum values 


e.e.. Ferrari et al.lll999a 

3). We 

. to Heeer et al.l (2003tl. 

SmarttI 


in detail. 

Numerically, we have /3ns = 0.024, fp = 0.136, $ns = 
5.97 X 10“®Mg^ and, consequently, Adns = 1.95 x 10~®Mg^. 

For the limits of integra tion ZfjuJ ) and zjvtf) 
in equation (ED: we follow iRegimbaul (1201 il l and 
iRegimbau fc Chauvineaul (l2007l l. namely 



if t'o < 


l+2;max 


— 1 otherwise; 


if Vo > 


l + Zn 


Zi{vo) = 


( 22 ) 


(23) 


— 1 otherwise, 


where l^min -— 10 Hz, i^max —^ 2 X 10 Hz, Ziain. — 0 
and Zmax = 20 are the minimum emitted frequency, the 
maximum emitted frequency, the minimum redshift and the 
maximum redshift, respectively. 


Using the ingredients shown above, integrating over 
the eccentricity e and substituting the numerical values, 
equation takes the form 


flgw(n,i/o) = 1.28 X 10"“ 


2/3 




„0.6{x-5.4) 


dz 


L 


0.333 -k [1 + zY/^E{z) 

5(n, £)/(£) 


F{e) 


-de, 


(24) 


where E{z) = iy0.3(l -I- zA + 0.7 and /(e) is given by either 
equation (HU or equation (O. 


6 RESULTS 

We show in Figs [T]l4] the spectra generated 
by some harmonics of the emitted signals, considering the 
distribution given by equation (nnj. It is worth mentioning 
that, although we show in these figures only some harmonics, 
we include in the calculations of the total spectra (Fig. 
all relevant harmonics. 

Note that for high eccentricities, the amplitudes 
corresponding to the higher harmonics become extremely 
relevant. Particularly, it is worth noticing the contrast 
between the cases corresponding to e = 0.25 and e = 0.9. 

We show in Fig. [5] the backgrounds generated when 
we sum the contributions of all relevant harmonics for the 
cases corresponding to Figs. lH4l We also include in this very 
figure the spectra for: (a) the uniform distribution given by 
equation (I17II : (b) DNS systems in circular orbits. 

Comparing the uniform curve in Fig. [5] with, for 
example, that by a population of DNSs in circular orbits, 
we can observe the effects of the orbital eccentricities: the 
frequency band of the background is broadened and shifted 
towards higher frequencies. 

The frequency band of the spectrum for circular systems 
ranges from ~ 10~® to ^ 10“®Hz, whereas for the uniform 
distribution function the spectrum ranges from ~ 10“® to 
~ IQ-'^Hz. 

Also, the amplitudes of the spectra are slightly reduced 
in the region between ~ 10“® and 10“®Hz whether 
the eccentricity of the population of DNSs is increased. 
This occurs due to the reduction of the amplitudes of 
the contributions corresponding to the lowest harmonics as 
compared to the highest ones. 

This effect can be more clearly observed in the other 
spectra present in Fig.O namely for e = 0.25, the maximum 
amplitude occurs at a frequency of ~ 10“®Hz; whereas for 
e = 0.9 the maximum occurs at ~ 10“'^Hz. As expected, 
the higher the eccentricity is, the higher the frequency 
corresponding to the maximum amplitude of the spectrum. 
Moreover, the amplitudes have maximum values of ~ 10“^® 
and ~ 10“^® for e = 0.25 and e = 0.9, respectively. 

A relevant issue is to investigate the detectability of 
the backgrounds studied here. Given the frequency bands 
involved, we consider the detectability by pulsar timing 
array (PTA) techniques and by the planned (space-based) 
antenna eLISA. _ 

We refer the reader to, for example. I Moore et al.l (l2015lf . 
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v(Hz) 


Figure 1. Backgrounds generated by some harmonics considering 
that all systems have eccentricity 0.25. 



v(Hz) 

Figure 2. Same as in Fig. 1 for eccentricity 0.5. 


where the sensitivity curves of several detectors, including 
the ones we mentioned here, are provided. 

Concerning the PTA technique, the frequency band 
of operation ranges from 10~^ to ~ 10~®Hz, with 
maximum sensitivity ~ 10~®Hz. Particularly, in terms of the 
energy density parameter, EPTA (IPTA) has a maximum 
sensitivity of ~ (~ 10“^^), whereas SKA is planned 

to have maximum sensitivity of ~ 

The planned eLISA (and LISA, its precursor) has 
maximum sensitivity, in terms of the energy density 
parameter, of 10”^°) around ~ 10“^Hz. 

Given the sensitivity curves of PTA and eLISA (LISA), 
we conclude that the backgrounds shown in Fig. [S] would not 
be detected. 


7 CONCLUSIONS 

In this paper we show how to obtain the stochastic 
background of GWs generated by long-lived DNSs in 
eccentric orbits. We calculate the spectra by means of a 
different method: we determine the spectra that would 
be generated by each harmonic and then we sum these 
contributions in order to obtain the total backgrounds. 



v(Hz) 

Figure 3. Same as in Fig. 1 for eccentricity 0.75. 



Figure 4. Same as in Fig. 1 for eccentricity 0.9. 


In addition, we study the role of the eccentricity on 
the resulting spectra. This is achieved by using different 
distribution functions for the eccentricity. Particularly, we 
use two distribution functions: the Dirac delta function given 
by equation (HU and the uniform distribution given by 
equation (HU. 

From Fig. [S] we note that the distribution of 
eccentricities affect the amplitudes and the frequency bands 
of the backgrounds. For example, we note that using the 
distribution given by equation 113, the amplitudes decrease 
in the frequency band ranging from ^ 10“® to 10“®Hz, 
relatively to the case of a population of systems in circular 
orbits (see Fig. o. Concerning the frequency bands of the 
backgrounds, they are enlarged and shifted towards higher 
frequencies. 

In general, we note that the higher the eccentricity of 
the population is the higher the frequency corresponding to 
the maximum amplitude of the spectrum. This behaviour 
is caused by the relative contribution given by the 
harmonics: for high eccentricities, higher harmonics affect 
more significantly the backgrounds. 

Concerning the detectability of the spectra studied here, 
they would not be detected by the presently proposed space- 
based antennas, since the amplitudes involved are very 
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Figure 5. Total spectra (obtained by summing all relevant 
harmonics) for the cases corresponding to Figs 1-4. Also shown, 
the backgrounds generated by DNSs in circular orbits and by the 
uniform distribution of the eccentricities given by equation II17II . 

small. Similar conclusion holds if we consider detections 
using the PTA technique. 

In forthcoming papers we will address the study to 
systems formed by two black holes and those formed by a 
black hole and a neutron star. In addition, we will consider 
the high-frequency regimes, in which we cannot neglect the 
effects of the time evolution of the orbital parameters. Such 
cases are worth studying because, among other things, the 
spectra generated could form foregrounds for the planned 
space based antennas eLISA, DECIGO and BBO. 

In addition, we expect that in the cases in which the 
time evolution of the parameters is considered, the form of 
the distribution f{e) will affect the results significantly, since 
in this case the frequency and the eccentricity will be related 
to each other. 
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